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We investigate the contribution of inelastic and elastic processes to single inclusive hadron pro- 
duction in proton-proton and proton (deuteron)-nucleus collisions at RHIC and the LHC. Using the 
hybrid formulation which includes both elastic and inelastic contributions, supplemented with the 
running-coupling Balitsky-Kovchegov equation, we get a good description of RHIC data. It is shown 
y—i • that inclusion of the inelastic terms makes the transverse momentum dependence of the production 

7—i ' cross section steeper in the mid-rapidity region but does not affect the cross section in the very 

forward region. The inelastic processes also lead to a sharper increase of the nuclear modification 
factor R p a with increasing pr- We also make predictions for the nuclear modification factor in 
q ' proton-nucleus collisions at the LHC (-^/s = 4.4 and 8.8 TeV) at various rapidities using the Color 

, Glass Condensate framework. 

Q: 

CN ' I. INTRODUCTION 

(N : 

i—i, The Color Glass Condensate (CGC) formalism [l| is a self-consistent, effective theory approach to QCD interactions 
• at high energy (or equivalently small x). Even though it is a weak coupling approach, it is different from the collinear 
Oh 1 factorization based approach of pQCD in two important aspects; first, it re-sums quantum corrections which are 
enhanced by large logarithms of 1/x as opposed to large logarithms of Q 2 in pQCD and second, it includes high gluon 
density effects which are important at small x and/or for large nuclei where the physics of gluon saturation may be 
. the dominant. 

The CGC formalism has successfully been applied to many QCD processes, from fully inclusive ones such as structure 
, functions in DIS to single and double inclusive particle production in proton-proton and proton-nucleus collisions at 
J> 1 high energy, see Ref. [l| and references therein. The CGC formalism has been also quite successful in providing 
predictions for particle multiplicities at the LHC 043] and may provide a first-principle way of understanding of 
isotropization and thermalization of QCD matter produced in high energy heavy ion collisions at RHIC and the LHC 

The observed suppression of single inclusive hadron production and the disappearance of the away side peak in 
double hadron production in the forward rapidity region of deuteron-nucleus (dA) collisions at RHIC are perhaps 
the strongest evidence for the importance and possibly dominance of saturation effects at RHIC. This will soon be 
further tested at the LHC where one will be able to probe CGC dynamics in a much larger kinematic region due to 
the larger energy of the collisions at the LHC. Single inclusive hadron production in proton-nucleus (pA) collisions 
at RHIC and the LHC has been investigated by many authors [8H13| in the CGC formalism using varying degrees 
of approximations and models (for an alternative description, see for example Refs. [3, EH). The most important 
ingredient of the single inclusive hadron production cross section which captures the saturation dynamics is the 
fundamental (or adjoint) dipole cross section, the imaginary part of the quark anti-quark scatte ring amplitude on a 
proton or nucleus target. This dipole cross section satisfies the JIMWLK/BK evolution equations [IsQjJ and re-sums 
the small x as well as high gluon density effects. The evolution equation for the dipole cross section is now known 
with next-to- leading-order (NLO) accuracy [l8j], see also Ref. flTjj . 

There are two distinct but related approaches to hadron production in high energy asymmetric (such as proton- 
nucleus or very forward proton-proton and nucleus-nucleus) collisions. One is the well-known kx factorized approach 
[2fl I2l1 ] where partons in both the projectile and target are assumed to be at very small x (x < 0.01) so that the 
CGC formalism is applicable to both the projectile and target. This approach is valid as long as one stays away from 
the projectile fragmentation region. An alternative approach was developed in [9| where one treats the projectile 
wave-function perturbatively, i. e. using the standard DGLAP picture while treating the target by CGC methods. 
This approach is better suited for the projectile fragmentation region. Very recently this approach has been improved 
by keeping the inelastic pieces of the cross section which may be important at high transverse momentum [221 ] . Here 
we numerically investigate the contribution of theses inelastic contributions to single inclusive hadron production in 
proton-nucleus and proton-proton collisions at RHIC and the LHC and show that inclusion of these terms improves 
the high px behavior of the cross section. The nuclear modification factor is also shown to increase faster with 
increasing px than the case where these inelastic contributions are ignored. 

The paper is organized as follows: In the next section, we introduce our formalism for the inclusive hadron production 
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in pA collisions, namely the hybrid formulation which includes both elastic and inelastic contributions, supplemented 
with the running-coupling Balitsky-Kovchegov equation. Section III is devoted to comparison with the experimental 
data and to a discussion of various predictions for the LHC energies. We conclude in Sec. IV. 



II. SINGLE INCLUSIVE HADRON PRODUCTION; MAIN FORMULATION 

The cross section for single inclusive hadron production in asymmetric collisions (scattering of a dilute system of 
parton with a dense one) at high energy is given by [22j 



dN pA^hX K 



dPpxdr) (27r) 
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where fj(x,Q 2 ) is the parton distribution functions (PDF) of the incoming proton which depends on the light-cone 
momentum fractions x and the hard scale Q. The function D h / i {z 1 Q) is the hadron fragmentation function (FF) of 
i'th parton to the final hadron h with a momentum fraction z. The variables 77 and px are the pseudo-rapidity and 
transverse momentum of the produced hadron. The longitudinal momentum fractions X\ and x 2 are defined as follows 

xf ~ — F e x\ = — ; x 2 = x\e (2) 
V s z 

where is the collision energy per nucleon. Here we neglect hadron masses since we are only interested in light 
hadron production at high-p-r (thereby rapidity and pseudo-rapidity are equal) . 

It is perhaps useful to remind the reader of the derivation of Eq. ([lj. The first line of Eq. (p} was first derived 
in Ref. @. The result of @ has been recently improved in [22j by keeping the inelastic pieces which lead to the 
second line in Eq. ([1]). Our main goal in this paper is to consider the effect of this new term in the inclusive hadron 
production at both RHIC and the LHC. Let us first focus on the first line of Eq. (p}, the DHJ term Q. We refer 
the reader to Ref. [9] for technical details and just outline the derivation of the DHJ term. One starts by calculating 
two particle production cross section in proton-nucleus scattering. The simplest process is when a quark from the 
projectile scatters on the target and radiates a gluon either before or after the scattering (see Fig. (12) in 23]). The 
incoming quark as well as the outgoing quark and the radiated gluon can all multiply scatter on the target. However, 
if one is interested in single inclusive production, one needs to integrate over one of the final state partons. Some of 
the Feynman diagrams will have collinear divergences whereas others do not. There is a collinear divergence in the 
final state which happens when the outgoing quark and the radiated gluon are nearly parallel and only the initial state 
quark multiply scatters on the target. This divergent term is absorbed into quark-hadron fragmentation function and 
lead to its evolution with Q 2 according to the LO DGLAP evolution equation. The finite (non-collinear divergent) 
pieces are ignored as they are part of the NLO corrections. 

There is also a collinear divergence in the initial state which happens when the incoming quark and the radiated 
gluon are nearly parallel and only the final state quark multiply scatters on the target. This sort of collinear divergence 
is absorbed into the incoming parton distribution function and leads to its evolution according to the standard DGLAP 
evolution equation. Again, the finite parts of these terms were ignored in Q since theycorrespond to higher order 
(in a s ) corrections. It was pointed out in p2j that the finite pieces which are ignored in Q may be important at high 
Pt and therefore can lead to a modification of the production cross section. Keeping the finite diagrams (which do 
not have a collinear divergence) and making the high px approximation (gradient expansion of the quadrupole cross 
section) leads to Eq. ([T]). While the above argument is for an incoming quark scatte ring on the target, inclusion of 
other processes such as an incoming gluon scattering on the target is straightforward [2J] where a similar analysis of 
the collinear divergences can be made. 

Now that the origin of these inelastic terms is made more clear, we comment on the relative significance of the two 
contributions. The first piece of eq. ([T]), dubbed the elastic part, corresponds to an incoming parton in the proton wave 
function scattering elastically on the target. This incoming parton initially has zero transverse momentum but picks 
up transverse momentum of order Q s after multiply scattering on the target. This term should be most important 
when the transverse momentum of the produced hadron is of order Q s or perhaps even a bit larger. The second term 
in eq. ([1} , dubbed the inelastic piece, corresponds to a high transverse momentum parton radiated from the incoming 
parton in the projectile wave function. This radiated parton is already at high transverse momentum and interacts 
with the target only once (higher number of scatterings will be power suppressed). This term is therefore important 
only when the produced hadron is at transverse momenta much higher that the saturation scale of the target Q s . 
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In Eq. (JTJ, we have introduced a K- factor to mimic the effect of higher order corrections. The inelastic weight 
functions Wj,/j are given by 
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where N c denotes the number of colors. The function P%u in Eq. (JTJ) denotes the Altarelli-Parisi splitting function 
that describes the probability of a given parton j splitting into two others. The leading-order splitting functions (for 
N c = 3) are given by j25j . 
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where rif is the number of active flavor and subscript + refers to the so-called "+" prescription used to regularize 
the singularities as £ — > 1 [25]. In Eq. (JXJ) , the amplitude Np (Na) is the two-dimensional Fourier transform of 
the imaginary part of the forward dipole-target scattering amplitude AfA(F) m the fundamental (F) or adjoint (A) 
representation, 



= I d 2 re- ikT - ? (l-M A{F) {r,Y = ln(zO/z))) 



(11) 



where r — \r\ is the dipole transverse size. In the large- N c limit, one has the following relation between the adjoint 
and fundamental dipoles, 



Na{t, y) — 2Af F (r, Y) — Afp(r, Y). 



(12) 



The amplitude Na(f) incorporates all multi-scatterings between a projectile color-dipole and the target and encodes 
the small-x d yna mics. In the CGC framework, it can be obtained from the solution of JIMWLK/BK evolution 
equations 0,113, an infinite set of coupled nonlinear equations for the different Wilson line correlators which system- 
atically incorporate small-x gluon emission to all orders In the large- 7V C limit, the JIMWLK evolution equations 
reduce to the Balitsky-Kovchegov (BK) equatio n 117 1, a closed-form equation for the evolution of the dipole amplitude. 

271 has the following simple form : 



The running coupling BK (rcBK) equation [Tt], l26l . 

dN A (F) (r, x) 



(Vh]! r , s - j d 2 ri K mn {r,r u r 2 ) [Af A ( F )(n, x) + N A (f){T2, x) - Af A (F){r, x) - Af A (F)(n, x) Af A (F) fa, x)] 
where the evolution kernel K run using Balitsky's prescription [27| for the running coupling is defined as, 
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(14) 



with fa = f—r\. For the running coupling in the above equation we use the scheme proposed in Ref. [34[ at one- loop 
level. Notice that in the master equations (|1I13[) . the impact-parameter dependence of the collisions was ignored. 
However, for the min-bias analysis considered here this may not be important. Nevertheless, it has been shown 
by several studies that impact-dependence of the BK equation is important at very large rapidities and for fixed 
centralities 0-H, H|l 29] , though it is very challenging to implement numerically. 
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III. NUMERICAL RESULTS AND PREDICTIONS 



We now evaluate the single inclusive hadron production cross section numerically. To do this, we will use the NLO 
MSTW 2008 PDFs [H and the NLO KKP FFs 31]. We have checked that AKK FFs H3 are also consistent with 
our results. We assume the factorization scale in the FFs and the PDFs to be Q = pt- 

The only input for the rcBK equation is the initial conditions for the evolution of the dipole amplitude which is 
commonly taken to be a McLerran-Venugopalan (MV) type model (33[: 



Af(r,Y=0) = 1 - exp 



In 



1 

Ar 



(15) 



where A = 0.241 GeV and 7 = 1.119 0, Notice that global fits to structure functions in DIS in the small-x 
region show that 7 > 1 is preferable and 7 ~ 1.119 provides a good fit to the DIS data [34). The onset of small-x 
evolution is taken at xq = 0.01. Then the only free parameter is the initial value of the saturation scale Qq s (probed 
by quarks) for proton and nucleus at xo — 0.01. The value of Qq s can be also fixed via a fit to the structure functions 
in electron-proton and electron-nucleus scatterings in the small-x region. In order to investigate the sensitivity of our 
results to the choice of the parameter sets, we take Qo s as a free parameter and let it be determined by the RHIC 
data in proton-proton (pp) and deuteron-gold (dAu) collisions. In this way, we will also examine whether a universal 
description of various low-x data can be achieved by the rcBK evolution equation. 

In order to facilitate a comparison of the solution to rcBK dipole evolution equation, we also employ the DHJ dipole 
parametrization [l(| for Naqa which has been used to describe the forward rapidity data in dAu collisions at RHIC. 
We refer the reader to Ref. [lC| for the details of the DHJ dipole model. The dipole scattering amplitude in the DHJ 
model is simply given by 
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where the anomalous dimension 7 DHJ in the DHJ model is parameterized as, 



7 UHJ (er,a;) = 7. + (l-7.) 
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with y — log(l/a;) and the scale Qt in the anomalous dimension is related to the inverse transverse size of the dipole 
s» 1/r 2 . Saturation scale Q a (x) is defined as Q 2 s {x) = A 1 Jj : 3 j : (x s q/x) x with A e // = 18.5 for the minimum bias 
dAu collisions. Here we are interested in investigating whether the new inelastic contribution, second line in Eq. ([1}, 
will affect the description of the RHIC data |10|. To this end, we take the same parameters for the DHJ dipole 
parametrization as employed in Ref. [l(| which provides a good description of the RHIC data without the presence of 
the inelastic contribution namely a s = 0. The parameters A = 0.3 and x R n = 10 -4 were extracted from a fit to HERA 
data, and parameter d was fitted to the RHIC data and set to d = 1.2 [Hj. The anomalous dimension in Eq. (fT?} 
runs from the LO BFKL value 7 S = 0.628 at small x to the DGLAP value 7 DHJ — > 1. This model incorporates the 
geometric scaling window as expected from the BK equation [35| consistently. 

Here our aim is not to fit the data but to use the best theoretical tools available in low-x physics to highlight 
the uncertainties involved in making robust predictions from the CGC formalism for the upcoming proton-nucleus 
collisions at the LHC. Therefore, we take K = 1 throughout this paper. 

In Fig. []] we show the single inclusive hadron production yields in pp and dAu collisions at RHIC y/s — 0.2 TeV 
at different rapidities using the DHJ parametrization of the dipole cross section as well as the rcBK dipole solution. 
In order to investigate the contribution of the inelastic term to single inclusive hadron production, we also show the 
results without this term, namely a s = in Eq. (p}, denoted as DHJ in Fig. [T] (left). The inelastic contribution term in 
Eq. |l]) is explicitly proportional to a s . Notice that in the derivation of this formula at the leading- twist order, a s was 
assumed to be a fixed parameter. It is not a priori obvious whether the running- coupling corrections to Eq. (|T|) can 
be simply incorporated by replacing a s to a running a s (Q). It has been shown for example that in fey-factorization 
formulation, the running coupling effect changes the equation (36|. We have checked that a s « 0.05 -7- 0.15 in Eq. ([T]) 
gives a reasonable description of RHIC data for both pp and dAu collisions. It is clear that inclusion of the inelastic 
terms improves the px dependence of the cross section closer to mid-rapidity while there is no visible contribution at 
the most forward rapidity considered. This is more clearly seen in Fig. [1] at the upper-left panel, where the inelastic 
contribution is seen to significantly improve the description of the data for more central collisions at 7/ = 1 and makes 
the p^-spectra steeper in agreement with the data. For more forward collisions, the available phase space is limited 
and inelastic contributions are consequently negligible independent of the value of the strong coupling a s . 
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FIG. 1: Single inclusive hadron production in proton-proton (upper panel) and deuteron-gold (lower panel) collisions at different 
pseudo-rapidities at RHIC obtained by the solution of the running-coupling BK equation, the so-called rcBK (right) and the 
DHJ (left) dipole model. Right: dashed and full lines refer to the results coming from the rcBK equation corresponding to two 
different initial values for the saturation scale at xo = 0.01. We have taken a s = 0.1 in Eq. {TJ for all curves. Left: dashed and 
full lines refers to the results when a s = (the DHJ term or the elastic contribution) and a s =0.1 (for the inelastic term), 
respectively. We have taken K — 1 in all panels. The experimental data are from Ref. [g]. 



As we have already pointed out, the value of Qo s at xq = 0.01 is the only free-parameter left to be fixed for a given 
solution of the rcBK equation. In Fig. [1] (right), we also show the effect of various choices for the initial saturation 
scale Qos- In the case of pp RHIC data, we found that values in Qq s = 0.168 -j- 0.336 GeV 2 range give a consistent 
description of data. However, a smaller value of Qq s = 0.168 GeV 2 may be more preferable, specially at very forward 
rapidities. This is understandable since the available phase space for multiple re-scattering is limited at very forward 
rapidity. Therefore, a lower initial saturation scale is required to describe the cross-section. We note that HERA 
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FIG. 2: Right: nuclear modification factor R p a for inclusive charged hadrons h + + h~ production at the LHC \fs = 4.4 TeV 
and r\ — 4 coming from the solutions of the rcBK with different initial values for the saturation scale (at xq = 0.01) for proton 
and nucleus. The dashed and full lines refer to the cases when the cross-section in both pp and pA collisions was obtained via 
Eq. (0 by taking a 3 = (only elastic contribution) and a 3 —0.1 respectively. Left: the scaled unintegrated gluon distribution 
Na(x,pt) x Pt as a function of transverse momentum pr at a fixed x = 10~ 5 obtained from the rcBK equation with two 
different initial values for the saturation scale Qos(xo = 0.01). 



data on proton structure functions prefer a lower value for the proton initial saturation scale 1 Qq s = 0.168 GeV 2 [Hj]. 
In the case of mini-bias dAu collisions, the initial nuclear (gold) saturation scale within Qq s — 0.5 -j- 0.67 GeV 2 is 
consistent with the RHIC data within the error bars. Unfortunately the available DIS data for nuclear targets are 
limited and have large experimental uncertainties. It is therefore difficult to pin down the exact value of Qq s for nuclei 
based on only DIS data, see also Ref. [37J. In the case of proton- nucleus collisions, due to theoretical uncertainties 
and rather large experimental data errors, it is also not possible to uniquely fix the initial value of Qq s - Nevertheless, 
the extracted values of initial nuclear saturation scale here are compatible with values extracted in other studies, see 
Refs. [3, [l2j and reference therein. 

It is seen from Fig.[T]that for the description of ir° production in both pp and dAu collisions at very forward rapidity 
(j) = 4) at RHIC, a iC-factor of ~ 0.4 ~ 0.6 may be needed. The necessity of such a small .ftT-factor at very forward 
rapidity at RHIC yfs = 200 GeV was also shown in Ref. [l2j where the inelastic contribution was ignored. Notice 
that as we have already pointed out, the effect of inelastic contribution at very forward rapidities at RHIC energy 
y/s = 200 GeV is negligible numerically, see Fig. Q] It is possible that at very forward ra pidi ties at RHIC energy, 
other mechanisms also partially contribute to the hadron production, see for example Refs. [38l [39|. 

Next, we present our predictions for single inclusive hadron production at the LHC in terms of the nuclear modifi- 
cation factor R p a hoping that some of the theoretical uncertainties, such as sensitivity to K factors, will be reduced. 
The Nuclear modification factor R p a is defined as 

R _ 1 dNP A ^ hx ^dNPP^ hx 
pA Ncoii cPprd-q cPprdri 

where N co u is the number of binary proton-nucleus collisions. We take N co u = 6.5, 7.4 at = 4.4 and 8.8 TeV, 
respectively 2 [40j |. 



1 The parameter set Qg 3 = 0.168 GeV 2 and 7 = 1.119 that we used in Fig. [T] also gives excellent description of the structure function 
data in e+p collisions with x 2 /d.o.f. = 1.104 l'34l . 

2 In order to compare our R p a predictions with experimental data, one may need to rescale R p a by matching the normalization N co u to 
the experimental value. 
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FIG. 3: Nuclear modification factor R p a for h + + h~ production in proton-nucleus collisions at the LHC (\/s = 4.4 TeV) at 
different rapidities (from top to bottom: r\ — 4, 6, 7) obtained from the solution of rcBK equation assuming two different initial 
nuclear saturation scales Qq s — 0.5 GeV 2 (left) and Qq s = 0.67 GeV 2 (right). In both panels the initial saturation scale for 
proton was taken Qq s = 0.168 GeV 2 . The effect of different value for the strong coupling a s in Eq. {l} is also shown. 



In Fig. [2] (right), we show the nuclear modification factor R p a for inclusive charged hadrons h + + h~ production at 
■s/s = 4.4 TeV and r\ — 4 obtained from different solutions of the rcBK equation corresponding to different values of 
Qos(^o = 0.01) extracted from RHIC data (see description of Fig. [1]). We also show the contribution of the inelastic 
term by showing the results due to only the DHJ term (a s = 0). The value of the strong-coupling in the inelastic 
term in Eq. (fTJ) is set to a s =0.1 (the same value was taken in Fig. [1]). It is obvious that taking different values 
for the saturation scale Qo s (xq = 0.01) for proton and nuclear targets significantly changes the nuclear modification 
factor. Therefore, the measurement of R p a provides vital information about the initial saturation scale of target and 
small-x evolution dynamics. Inclusion of the inelastic term changes R p a and makes it increase faster at high-pT, see 
also Fig. |3] Notice that rcBK solutions taken here approximately reproduce the perturbative power-law behavior of 
the dipole-amplitude Na(f) ~ V.Pt a * high-PT, see Fig. [5] (left). We recall that the parameters of rcBK solutions 
used here were obtained from a fit to HERA data for virtuality Q 2 E [0.25,45] GeV 2 [H. Therefore, our results at 
very high-py may be less reliable. 

In Fig. 21 we show our predictions for R p a for h + + h~ production at ^fs = 4.4 TeV and 77 = 4 — 7 using the solution 
of the rcBK dipole evolution equation, Eq. (TIB"]) , assuming initial nuclear saturation scales of Qq s — 0.67 GeV 2 (right 
panel) and Qq s — 0.5 GeV 2 (left panel). In both panels we have assumed the initial saturation scale of proton to 
be Qq s — 0.168 GeV 2 . We note that a larger initial saturation scale for the nucleus leads to a faster rise of R p a 
with transverse momentum . For comparison, in Fig. |4] (right panel) we show the corresponding R p a obtained using 
the DHJ dipole model, defined in Eq. (jTHJ) . It is seen that both approaches lead to a suppression of R p a at forward 
rapidities at the LHC and that the DHJ parameterization leads to a flatter transverse momentum dependence. We 
recall that both the rcBK solution and the DHJ model provide a reasonable description of RHIC data. 

In order to highlight the uncertainties associated with the different choices of the strong-coupling constant in Eq. (Q]) 
more clearly, in Fig. 2] (left) we show R p a for three different values of a s namely, a s — corresponding to the elastic 
term only, and a s — 0.1,0.15 for the inelastic contribution. It is seen from Figs. [2 [3l [4] that a t rapidities close to 
mid-rapidity, increasing a s reduces R p a while at very forward rapidities and high-p-^ the opposite happens. 

In Fig. [S] (right), we show our predictions for R p a for inclusive charged hadron production at yfs — 8.8 TeV and at 
different rapidities obtained from the rcBK equation (|13[) with different values of the strong coupling in the master 
equation ([!]). It is seen that the energy-dependence of R p a from y^s = 4.4 to 8.8 TeV is rather weak 3 . From Figs. El 



3 Note that our results for R p a a t the LHC without inclusion of the inelastic term is different from Ref. [T^l mainly due to the fact that 
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P T [GeV] p T [GeV] 

FIG. 4: Right: R p a for inclusive charged hadron production at different rapidities at the LHC obtained from the DHJ 
parameterization of the dipole profile Eq. (|16[) . with different values for the strong coupling in Eq. ((T]). R p a for inclusive 
charged hadron production for various values of the strong coupling constant a a in Eq. (JTJ at the LHC (y/s = 4.4 TeV and 
r\ = 5) obtained by the rcBK equation Eq. (|13p . 



21 we also note that at very forward rapidities the uncertainty associated with the choice of a s is reduced. This is 
in accordance with the fact that the effect of inelastic contribution at very forward rapidities is negligible, see also 

Fi g . m 

Our prediction for 77 = at the LHC energy 4.4 TeV is shown in Fig. [3] (left), using the solution to the rcBK 
evolution equation and assuming two different initial nuclear saturation scales of Qq s = 0.5, 0.67 GeV 2 (extracted 
from RHIC data). In Fig. [SJ we assumed the initial saturation scale for proton to be QfL = 0.168 GeV 2 (extracted 
from RHIC and HERA data) . The theoretical error bars in Fig. [5] show the uncertainties mainly associated with 
the choice of a s in Eq. |T]). The observed suppression of R p a at midrapidity and high px for the case of lower 
initial nuclear saturation scale QcL = 0-5 GeV 2 is larger compared to the results obtained within the gluon saturation 
approach with quasi-classical (Glauber) approximation [l3j]. We note that there are large uncertainties in R p a in 
midrapidity at the LHC due to the choice of the initial saturation scale for the rcBK evolution equation (fl~3f , and the 
value of strong-coupling constant in Eq. (TTJ. More importantly, large sensitivity of R p a to the value of a s in Eq. (TTJ) 
in midrapidity at the LHC indicates that higher order corrections in Eq. |T]) should be important in midrapidity at 
the LHC energy. Therefore, we believe that our predictions for R p a at midrapidity may be less reliable compared to 
our results for the very forward rapidity collisions. 

It should be noted that the particle production cross-section given by Eq. |T]) is intrinsically asymmetric, namely it 
treats the projectile proton approximately in the collinear factorization framework while treating the target proton (or 
nucleus) in the CGC framework. Strictly speaking, this may be justified only for particle production in the collision 
of a dilute system on a dense system, such as particle production in mid or forward rapidity in pA collisions or 
in particle production in the very forward rapidity region in symmetric collisions, such as proton-proton or nucleus- 
nucleus collisions. Therefore, our formalism can not be reliable for particle production in midrapidity in proton-proton 
collisions. Unfortunately, this is also what one needs in order to calculate the nuclear modification factor R p a in 
midrapidity. A better approach to particle production in midrapidity in symmetric collisions where both projectile 
and target are dilute (for p t » Q s ) might be to use the k t factorization formalism, proven to LO accuracy for single 



inclusive hadron production 41]. This is beyond the scope of the present work and we leave it for a future study. 



we have used a different parameter set for the rcBK solution with 7 = 1.119 which was recently suggested in Refs. 13, |34H . For the 
sensitivity of R p a to the various allowed solutions of the rcBK equation, see Fig. [2] and related discussions. 
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FIG. 5: Right: nuclear modification factor R p a for h + + YT production in proton-nucleus collisions at the LHC {^fa = 8.8 
TeV) at different rapidities obtained from the solution of rcBK dipole evolution equation (|13|l with different values for the 
strong coupling constant in Eq. (JTJ. Left: R p a for h + + h~ production in proton-nucleus collisions at the LHC in midrapidity 
(y/s = 4.4 TeV, n — 0) for two different initial nuclear saturation scales of Qq s = 0.5, 0.67 GeV 2 extracted from RHIC data. The 
initial saturation scale for proton is taken to be Qq s = 0.168 GeV 2 . The theoretical error bars mainly show the uncertainties 
associated with the choice of a s . 



IV. SUMMARY 



We have quantitatively studied, for the first time, the contribution of both elastic and inelastic processes to single 
inclusive hadron production cross section at RHIC and the LHC using the CGC formalism. We observe that inelastic 
contributions to single inclusive hadron production are significant at high transverse momentum and close to mid- 
rapidity. On the other hand, their contribution is very small in the forward rapidity region. Furthermore, we note that 
inclusion of these inelastic terms makes the nuclear modification factor R p a grow faster with increasing transverse 
momentum. We make detailed predictions for R p a at the LHC using the numerical solution of the running-coupling 
BK equation. We have studied various theoretical uncertainties associated with the choice of the initial saturation 
scale Qq s for a proton and nucleus. We have shown that the nuclear modification factor R p a measured at the LHC in 
the forward rapidity region is a sensitive probe of the low-x dynamics and can help constrain Qq s further. We have 
shown that various theoretical uncertainties in our formalism are minimized at very forward rapidities at the LHC. 
Therefore, measuring the nuclear modification factor R p a m the very forward region in proton-nucleus collisions at 
the LHC will be a robust test of gluon saturation dynamics and the Color Glass Condensate formalism. 
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